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Abstract 


Much of the research on the behavior of the SIS model on networks has 
concerned the infinite size limit; in particular the phase transition between 
a state where outbreaks can reach a finite fraction of the population, 
and a state where only a finite number would be infected. For finite 
networks, there is also a dynamic transition—the immortality transition— 
when the per-contact transmission probability A reaches one. If A < 1, the 
probability that an outbreak will survive by an observation time t tends to 
zero as t —> oo; if A = 1, this probability is one. We show that treating A = 
1 as a critical point predicts the A-dependence of the survival probability 
also for more moderate A-values. The exponent, however, depends on 
the underlying network. This fact could, by measuring how a vertex’ 
deletion changes the exponent, be used to evaluate the role of a vertex 
in the outbreak. Our work also confirms an extremely clear separation 
between the early die-off (from the outbreak failing to take hold in the 
population) and the later extinctions (corresponding to rare stochastic 
events of several consecutive transmission events failing to occur). 
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1 Introduction 

The Susceptible-Infectious-Susceptible model captures the dynamics of an in¬ 
fectious disease spreading in a population where infected people, upon recovery, 
become susceptible again. It has a long history, both as practical tool for pre¬ 
dicting and understanding real outbreaks (e.g. Ref. ID). and as a problem in 
applied mathematics (as a special case of stochastic logistic processes) [21 [3, [U 
[misnziisnm ini ED Cl] , statistical physics and computational science [51 [HI fTUl 
\mmm- The theoretical development has, for example, focused on fully con¬ 
nected topologies (or well-mixed models in epidemiological jargon) where every 
individual is equally likely to meet everyone else, at every unit of time. One 
of the recent advances is to calculate the exact value of the average extinction 
time in a finite system mm]- In parallel, one of the greatest advances of 
theoretical epidemiology the last two decades is to move away from the well- 
mixed assumption and model the contact over which the disease spreads as a 
network PU0Q3]. In this direction, researchers have, for example, calculated 
the epidemic threshold for a given network ® H CUE! — i.e. the critical value of 
the per-contact transmission probability A, below which a disease cannot reach 
a finite fraction of the population in the TV —> oo limit. 

For finite sized networks, there is also a phase transition, at least a thresh¬ 
old phenomenon albeit one whose value is trivial. If A < 1, any outbreak will 
eventually die out in a finite network; if A = 1, they will live forever. This im¬ 
mortality transition may be an epidemiological curiosity, as assumptions, such 
that the underlying network is relatively static, would break down for large 
enough A-values. It is also trivial in the sense that it is not an emergent phe¬ 
nomenon. Still, it could be the case that the increasingly unlikely extinction 
events as A grows follow the same statistics as the critical slowing down around 
thermodynamic phase transition. In other words, we cannot a priori exclude 
the possibility of something like a critical behavior of the immortality transition. 
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Figure 1: The two networks that we study (with the tree vertices of largest 
devalues marked). Panel (a) shows a sexual network from the beginning of the 
American HIV outbreak of the late 1970’s from Ref. [19] ; (b) shows a sexual 
network of Icelandic MSM from Ref. EH. 


So, since we cannot rule it out, we will assume that it is indeed true and see 
where it leads us. 

In the remainder of this paper, we will investigate the scaling of the survival 
probability £—the chance, as a function of A, an outbreak survives past an ob¬ 
servation time t. Furthermore, we will investigate the contribution of individual 
vertices to the behavior of £. This, we will argue, gives a new way of looking at 
importance of vertices in the SIS model on networks. 


2 Preliminaries 

We use the constant-duration version of the SIS model defined as follows. As¬ 
sume an underlying network represented as a simple graph G = (V,E), where 
V is a set of N vertices and I? is a set of M edges. At the beginning of the 
simulation, all vertices are susceptible except one vertex that is infectious (or 
rather, becomes infectious at this very time step). Then, for every edge be¬ 
tween an infectious and a susceptible, the susceptible will become infectious 
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with a probability A the following time step. An infectious vertex becomes 
susceptible again 5 time steps after becoming infectious. When there are no 
infectious vertices, the outbreak is dead. For large enough, 6 this is effectively 
a one-parameter model -the product XS determines the entire behavior of the 
model, their actual values only set the time scale. If S is smaller than the short¬ 
est time for an outbreak to spread through the graph, there might be other 
effects arising. Ref. [TfT] discusses this issue and argues that for practical pur¬ 
poses the constant-duration version of compartmental models is equivalent to 
the constant-recovery rate version (the latter being more common in the mathe¬ 
matical literature). (But a disclaimer is that Ref. [U>j did not deal with extreme 
events such as extinctions.) We use 5 = 5, which is around the radius of the 
graphs we study (and thus of the order of the fastest times of the disease to 
spread across the network). 

Our core quantity for monitoring the outbreaks is the survival probability 
/(A, f)— the fraction of outbreaks that are dead by time t. df/dt gives prob¬ 
ability distribution of extinction times, which is perhaps a more a commonly 
studied property |3 El 0 Eld. All simulations are averaged over 10 6 indepen¬ 
dent simulation runs. 

For the purpose of this paper, any kind of small and somewhat heteroge¬ 
neous network would suffice. But since there are such empirical networks in the 
epidemiology literature, we will take two of them as our study objects. Both 
networks represent sexual networks of men who have sex with men (MSM) and 
were collected to study the HIV epidemics. Of course, other compartmental 
models are more appropriate than SIS for modeling HIV but there are 
other infections that spread over these networks— like chlamydia, syphilis and 
gonorrhea—that fits the SIS picture. Since these infections have different trans¬ 
mission pathways among MSM, epidemiologists typically treat MSM as a special 
case; also for diseases other than HIV m- 

The first example network is based on contact tracing (following infection 
chains backward in time) in the early HIV outbreak in Canada and the United 
States |T2] (we refer to this data as America). It has N = 40 and M = 41. 
Note that contact tracing induces structures in the sampled networks that do 
not have to be present in the underlying sexual networks [2D]- Then again, 
this networks serves well as an example. Another network represent the MSM 
network of Iceland in the early 1990 ’s and probably suffers less from sampling 
bias than the previous network m- Here, N = 75 and M = 115. We refer to 
this data Iceland. Both networks are visualized in Fig. [T] 


3 Results 

3.1 The scaling of the survival probability 

We start our expose of numerical results by plotting the survival probability 
for fixed observation times t as a function of A (Fig. [2]). This type of plot is 
reminiscent of the average outbreak size 9 as a function of A, which for infinite 
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Figure 2: In panels (a) (for the America data) and (b) (for the Iceland data), 
we see the survival probability £ as a function of A for an exponentially growing 
set of observation times t = 200 x 2 V in (a) and t = 50 x 2" in (b), where 
v = 0,..., 9 for both panels. The lower t value a curve has, the more to the left 
it is. 
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Figure 3: Fits of the extinction time data to an exponential functional form to 
determine the time constant r. Panel (a) shows graphs for America ; panel (b) 
is a corresponding plot for the Iceland. Error bars are smaller than the symbol 
sizes. 


6 





A A 

Figure 4: The time constant of the decay of £, t, as a function of A. Panel (a) 
gives results for the America data; panel (b) for the Iceland data. The curves 
are fits to the form Aexp(A/l) + B( 1 — A) - ^. The parameter values for the fit 
are A = 0.94(7), A = 0.0307(8), B = 2.9(5) x 1CT 4 and < = 78.6(6) for America , 
and A = 1.22(2), A = 0.0132(1), B = 8.2(4) x 1CT 4 and C = 197(6) for Iceland 
(the digits in the parentheses are the standard errors in order of the last digits). 


systems show a threshold phenomenon, f l (precisely defined as the fraction of 
individuals that at some time are infected) changes, at a critical A-value, from 
zero to > 0. At a first glance, it seems like the same thing happens in Fig. [2] 
At a certain A value, £ increases very rapidly from a value close to zero. The 
increase becomes steeper the larger the i-value is, so t seems to take the role of 
N in finite-size analysis of the phase transition in f2. As mentioned, we already 
know that there is a threshold behavior in this case, but at A = 1. At a closer 
look, we can see that the rise of the ^-curves happens later for every larger 
t-value. Since we use an exponential progression of t-values this increase is in 
effect very slow. Another conspicuous feature of these curves is that they all 
group into one for large enough f-values. This envelope of curves (or, rather, 
pseudo envelope, since it will disappear as t —>• oo) corresponds to the early die¬ 
off seen in many epidemic models. More precisely, one minus the £ value of the 
envelope gives the fraction of outbreaks that fail to take hold in the population. 
The remaining runs do reach some quasi-endemic state, but will eventually die 
due to the fluctuations in a finite network. 

Next we look closer at £ as a function of t for fixed A-values. See Fig. [3} We 
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chose A-values to get a large range in £ (i.e. around where curves are the steepest 
in Fig. [2]), but for both networks and all A the observation is the same—except 
very small t, £(f) follows an exponential function exp(— t/r) quite accurately. 
This has been found analytically for fully connected networks Hum, and we 
guess that it holds for all connected networks. The exact functional form of the 
early die-off (where £(f) decays faster than an exponential) could depend both 
on the heterogeneities of the networks around the seed and the discrete time. 
We will not go deeper into this, but focus on r. There are plenty of articles 
on the behavior of the average r in fully connected networks in the large -N 
limit MM- Perhaps the theories in these references could be extended to 
the finite-TV and derive r(A). As mentioned, we take a computational physics 
approach. Our next step is to construct a scaling ansatz for r(A). Assuming two 
independent scaling regimes where the large-A one is dominated by fluctuations 
similar to a critical point from the A = 1 transition—r ~ (1 — A) - ^ where 
£ is an exponent corresponding to a critical exponent. For the low-A part, 
the fully connected, N —> 00 case predicts a linear dependence [3], but that 
does not fit our numerical results, which seem to be exponential. We combine 
(more precisely, add, assuming they stem from independent mechanisms) this 
observation with the finite-size scaling form to get 

r = Aexp(A/Z) + B(1 — A) - ^ (1) 

where A and B are constants, l controls the small A behavior and £ determines 
the dynamics close to the immortality transition. Fitting to the form given by 
Eq. 0 is very accurate (Fig. |4j) . Fits to other four-parameter functional forms 
(with a linear or power-law scaling of the low-A term), is visibly worse (failing 
at low A). Unfortunately, only limited regions of A are accessible—for small A, 
the disease dies too fast to get reliable data; for large A, the simulations take 
too much time. The accuracy of the scaling form is probably higher than many 
numerical phase transition studies, which makes this somewhat more tolerable. 
Our conclusion from Fig. [4] is that it is consistent with a critical scaling around 
the A = 1 transition affecting the dynamics at much lower A. At the same 
time, we hope future studies could derive both terms of Eq. (|T|) in a more 
systematic way. As a final note, the values of £ are both high (compared to 
critical exponents) and different for the two networks. These both observations 
suggest that the network topology determines the exponent. This is different 
from the universality classes of spin systems on lattices where the exponents 
take rather few values. 

3.2 Vertex-vise contribution to £ 

We have seen that £ depends on the network, but how? Now we turn to examin¬ 
ing the contribution of each vertex to the scaling behavior of £. We quantify the 
contribution of vertex i by the ratio a between £ after and before i is removed. 
Let Gi represent G with i deleted, then: 


ch = £(Gi)/£(G) 


(2) 


Measure 

America 

Iceland 

a 

ki 

0.73(4) 

0.974(2) 

o3 

S-i 

rii 

0.82(4) 

0.75(5) 

a 

mi 

0.83(3) 

0.965(2) 

o 


0.64(4) 

0.917(6) 

a 

max Ki 

0.76(5) 

0.98(2) 


for a 

0.17(8) 

0.038(5) 

03 

a 

maxi?j 

0.72(6) 

0.97(4) 

T—1 

for d 

0.99(1) 

0.99(1) 


Table 1: Pearson’s correlation coefficient between a,;—capturing a vertex’ in¬ 
fluence on £—and standard descriptors of the position of a vertex in a network. 
The number in the parenthesis is the standard error (obtained through jackknife 
resampling) in the order of the last digit. For the one-parameter quantities we 
also show the optimizing parameter values (a for Katz centrality and d for 
PageRank). 


For all the vertices in our networks, £(G*) is less than £(G), or within one 
standard error from it. It is easy to understand that adding a vertex is like 
making the road to extinction one step longer and thus even harder to reach, 
but how this translates to the more dramatic divergence for large A is not clear. 
We leave the observation that a larger graph has larger ( as a conjecture and 
focus on the network structural predictors of at. 

First, we investigate parameter-free descriptors of V s position in the net¬ 
work, such as degree fcj, the number of vertices n, and edges mrii in the largest 
connected component of Gi, coreness, closeness, betweenness, and the eigenvec¬ 
tor centrality e^. The latter comes from the idea that a vertex’ centrality is, 
recursively, proportional to the sum of its neighbors’ centralities. This leads to 
ei being the i’th element of G’s adjacency matrix (where the element on the i’th 
row and j’th column is one of i and j form an edge, and zero otherwise). All the 
quantities we use are discussed in network theory textbooks such as Ref. [ T2\ . 
We will consider the ones that have the highest explanatory power (measured 
through the absolute value of Pearson’s correlation coefficient r) with respect to 
a,;, namely k t , rij. rrii. and e,;. The r -values between these measures and a are 
presented in Table [l] Even though the quantities that we do not list (closeness, 
etc.) are consistently worse, the two data sets show quite different results for 
the listed quantities. Note, from Fig. [I] that the Iceland network is still con¬ 
nected into one large component no matter which vertex that is deleted. In the 
America data, on the other hand, deleting the vertices with the largest cq values 
disconnects the network. For this latter data set, the quantities that capture 
the fragmentation of the network, i.e. rq and irii, are the best predictors of a^. 
In Fig. |T] we plot the three vertices of largest cq for the two networks, which 
illustrates this point well. For the other network, Iceland , degree and mi (that 
in this case—with little fragmentation), is strongly correlated with degree, are 
the ones with highest correlation with cq. Interestingly eigenvector centrality 
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is performing poorly, highlighting the difference between indirect interaction in 
the SIS model and the linear coupling that the eigenvector centrality builds on. 
We note the ranges of ai is [1.06(7), 1.9(1)] for America and [1.004(2), 1.24(5)] 
for Iceland (where the numbers in parentheses are standard errors in the order 
of the last digit). 

In addition to the zero-parameter importance measures above, we also mea¬ 
sure the best possible correlation for two one-parameter measures—Katz cen¬ 
trality Ki and PageRank I?.;. These measures are related to the eigenvector 
centrality—Katz centrality also assume the centrality is proportional to the sum 
of the centrality of the neighbors, but also plus a constant a for the vertex itself. 
PageRank is proportional to the occupation probability of an unbiased random 
walker that with a probability d moves to a neighbor of its current vertex, and 
otherwise it moves to a random vertex. The fact that a disease can spread from 
one vertex to many others make Katz centrality seem more appropriate. We in¬ 
clude PageRank for comparison. The optimal correlation with Katz centralities 
is indeed always stronger than for the PageRank. The PageRank is optimized 
in the small-d limit and the Katz centrality for an intermediate a value. For the 
Iceland data the Katz centrality shows the highest correlation of all, meaning 
that there are measurable higher order structures that captures at better than 
degree. For the America data, rij and to,; are still the measures with highest 
correlation with but after those come Katz centrality. In summary, for our 
somewhat sketchy analysis, how much deleting a vertex would fragment the 
network seems to be the most important structure for explaining a,;, degree is 
the second most important (but not the only other) factor. 

4 Discussion 

We have numerically investigated extinction events in the SIS model on small 
networks. Our observations are consistent with the assumption that the transi¬ 
tion at A = 1—below which an outbreak would always die out in a finite-sized 
network—can be treated with standard finite-size scaling theory, but with the 
size replaced by time (not an entirely new idea, cf. Ref. [13]). Using this as¬ 
sumption, we find exponents that are dependent on the network topology. The 
fact that the £ depend on the topology is different from critical phenomena 
where exponents are groped into universality classes. One interpretation is to 
see the network as an integral part of the model—after all, we we have already 
accepted to drop size-scaling from the picture. On the other hand, we could see 
this as an indication not to push the analogy between the immortality transition 
and critical behavior further. We argue that £ can be used as a parameter-free 
index of a vertex’ role in the SIS extinction dynamics. We define the index as 
the exponent for the network without the vertex, divided by the exponent with 
the vertex present. This index, we show, depends much (but not only) on how 
deleting the vertex would fragment the network. 

The extinction time for large A is extremely long. Even though our networks 
are small, if A > 1/2, one would have to wait more than 10 10 ) time steps for 


10 


the survival fraction to go below 50%, even in the smallest networks. In our 
analysis, a time step represents a fifth of the duration of the disease. This 
means that for diseases lasting a week or so, we would have to wait 10 5 times 
the age of the universe (and many more times an average lifespan, let alone 
the typical lifetime of an edge in the network) to see more than half of the 
outbreak die out. In addition, the extinction times grow fast with N (for fully 
connected networks, the growth of the average is exponential 0 )- From this 
discussion, we understand that no real outbreak would ever come close to the 
immortality transition without violating the assumption of a stable underlying 
network. Indeed, no direct SIS simulation comes close either. 

For practical purposes, the most important observation is the extremely 
clear separation between the early die-off of outbreaks that fail to take hold 
in the population, and the later extinction events from rare stochastic events. 
This is maybe most clearly visible as the envelope of the curves in Fig. [f] This 
observation has been made before—e.g. Ref. [3j argues for the relevance of the 
quasi-stationary stage of the SIS dynamics. Since awareness and countermea¬ 
sures can, effectively, rewire the network fairly quickly only the early extinction 
events are of practical interest. The fraction of early extinctions is thus, from 
a modeling point of view, very well defined, which makes it easier to make 
conclusive statements from simulations. From a theoretical point of view, the 
crossover from exponential to power-law divergence scaling (Eq. 0 ) is also 
interesting and opens for future investigations. 


Acknowledgments 

The author acknowledges suggestions from Beorn Jun Kim, Mark Newman, Jari 
Saramaki, and Taro Takaguchi. The author was supported by Basic Science 
Research Program through the National Research Foundation of Korea (NRF) 
funded by the Ministry of Education (2013R1A1A2011947). 


References 

[1] H. W. Hethcote and J. A. Yorke, Gonorrhea transmission dynamics and 
control (Springer, Berlin, 1984). 

[2] R. H. Norden, Adv. Appl. Prob. 14, 687 (1982). 

[3] R. J. Kryscio and C Lefevre, J. Appl. Prob. 26, 685 (1989). 

[4] I. Nasell, J. Theor. Biol ., 211, 11 (2001). 

[5] C. R. Doering, K. V. Sargsyan, and L. M. Sander, Multiscale Model. Simul. 
3, 283 (2005). 

[6] C. R. Doering, K. V. Sargsyan, L. M. Sander, and E. Vanden-Eijnden, J. 
Phys.: Condens. Matter 19, 065145 (2007). 


11 


[7] O. Ovaskainen and B. Meerson, Trends Ecol. Evol. 25, 643 (2010). 

[8] P. van Mieghem, EPL 97, 48004 (2012). 

[9] C. Castellano and R. Pastor-Satorras, Phys. Rev. Lett. 105, 218701 (2010). 

[10] G. Brightwell and M. Luczak, e-print arXiv:1312.7449 (2013). 

[11] M. Boguna, C. Castellano, and R. Pastor-Satorras, Phys. Rev. Lett. Ill, 
068701 (2013). 

[12] R. Pastor-Satorras, C. Castellano, P. Van Mieghem, and A. Vespignani, 
e-print arXiv: 1408.2701 (2014). 

[13] S. Gomez, A. Arenas, J. Borge-Holthoefer, S. Meloni, and Y. Moreno, EPL 
89, 38009 (2010). 

[14] M. Morris, Network epidemiology: A handbook for survey design and data 
collection (Oxford University Press, Oxford, 2004). 

[15] M. J. Keeling and K. T. D. Eames, J. R. Soc. Interface 2, 295 (2005). 

[16] P. Holme, Journal of Logistical Engineering University 30, 1 (2014). 

[17] J. S. Koopman, et al., J. Acquir. Immune Defic. Syndr. 14, 249 (1997). 

[18] P. D. Benn, et al., Sex. Transm. Infect. 83, 106 (2007). 

[19] D. M. Auerbach, W. W. Darrow, H. W. Jaffe, and J. W. Curran, Am. J. 
Med. 76, 487 (1984). 

[20] K. T. D. Eames and M. J. Keeling, Proc. R. Soc. Lond. B 270, 2565 (2003). 

[21] S. Haraldsdottir, S. Gupta, R. M. Anderson, J. Acquir. Immune Defic. 
Syndr. 5, 374 (1992). 

[22] M. E. J. Newman, Networks: An Introduction (Oxford University Press, 
Oxford, 2010). 

[23] M. J. Lee, S. D. Yi, and B. J. Kim, Phys. Rev. Lett. 112, 074102 (2014). 


12 




